In [1]:
% matplotlib inline
import numpy as np
from PyQt4 import QtGui
import sys
import os
import matplotlib.pyplot as plt
import mpl_toolkits.basemap.pyproj as pyproj
from pydelft import grd, dep
import datetime
import math

load a grid


In [2]:
grid = grd.grd()
grid.read_grd(r'/Users/spmls/GitHub/pydelft/dev/test_data/test.grd')
print('m: %s\nn: %s' % (grid.m, grid.n))


m: 11
n: 11

load a depth file


In [3]:
depth = dep.dep()
depth.read_dep(r'/Users/spmls/GitHub/pydelft/dev/test_data/test.dep')
print('depth shape: (%s, %s)' % (np.shape(depth.depth)[0], np.shape(depth.depth)[1]))


depth shape: (11, 11)

plot


In [4]:
plt.pcolormesh(grid.x, grid.y, depth.depth)


Out[4]:
<matplotlib.collections.QuadMesh at 0x10c9374a8>

create a depth file that plots on the grid shown above


In [5]:
Z = np.random.rand(grid.m, grid.n)

Z = np.array([np.append(i,-999.) for i in Z])
Z = np.insert(Z,-1,np.ones(np.shape(Z[0]))*-999.)
Z = Z.reshape((grid.m+1, grid.n+1))

np.savetxt(r'/Users/spmls/GitHub/pydelft/dev/test_data/random_test.dep', Z, delimiter = '\t', fmt = '%.3f')

load the depth file that was created


In [6]:
depth2 = dep.dep()
depth2.read_dep(r'/Users/spmls/GitHub/pydelft/dev/test_data/random_test.dep')

plt.pcolormesh(grid.x,grid.y,depth2.depth)


Out[6]:
<matplotlib.collections.QuadMesh at 0x10ca700f0>

test the pydelft function


In [7]:
Z = np.random.rand(grid.m, grid.n)
depth3 = dep.dep()
depth3.write_dep(Z, fname = r'/Users/spmls/GitHub/pydelft/dev/test_data/test2.dep', grid_fname = grid.filename)
plt.pcolormesh(grid.x, grid.y, depth3.depth)


saved depth file: test2.dep
Out[7]:
<matplotlib.collections.QuadMesh at 0x10cacd9e8>

test without a grid already made


In [8]:
Z = np.random.rand(15,4)
depth4 = dep.dep()
depth4.write_dep(Z, fname = r'/Users/spmls/GitHub/pydelft/dev/test_data/test4.dep')


saved depth file: test4.dep

make a grid for the above depth file


In [9]:
coord_system = 'Cartesian'
m = 15
n = 4
cellsize = 10
x0 = 0
y0 = 0

records = []
rows = math.floor(n/5)
remainder = n%5

records.append('* \n* Delft3d- rectilinear grid file created with pydelft \n* File creation date: %s\n* \n' \
%str(datetime.datetime.now()))
records.append('Coordinate System = %s\n' % coord_system)
records.append('\t%i\t%i\n' % (n,m))
records.append('0 0 0\n')

# Values for x
for i in range(1,m+1):
    etax = np.arange(x0, n*cellsize, cellsize)

    if n > 5:
        etax_1 = etax[0:5]
        records.append('ETA=\t%i\t%.20e\t%.20e\t%.20e\t%.20e\t%.20e\n' %(i, 
                                                                         etax_1[0], 
                                                                         etax_1[1], 
                                                                         etax_1[2], 
                                                                         etax_1[3],
                                                                         etax_1[4]))
        etax_mid = etax[5:-remainder]
        for k in np.arange(0,np.size(etax_mid), 5):
            records.append('      \t\t%.20e\t%.20e\t%.20e\t%.20e\t%.20e\n' % (etax_mid[k],
                                                                              etax_mid[k+1],
                                                                              etax_mid[k+2],
                                                                              etax_mid[k+3],
                                                                              etax_mid[k+4]))
        etax_last = etax[-remainder:]
        records.append('      \t\t%.20e' % etax_last[0])
        for i in etax_last[1:]:
            records.append('\t%.20e' %i)
        records.append('\n')
    elif n <= 5:
        etax_1 = etax
        records.append('ETA=\t%i' % i)
        for i in etax_1:
            records.append('\t%.20e' % i)
        records.append('\n')

# Values for y
y = np.arange(y0, m*cellsize, cellsize)
for j in range(1,m+1):
    if n > 5:
        etay = np.ones((n,1))*y[j-1]
        etay_1 = etay[0:5]
        records.append('ETA=\t%i\t%.20e\t%.20e\t%.20e\t%.20e\t%.20e\n' %(j, 
                                                                    etay_1[0], etay_1[1],
                                                                    etay_1[2], etay_1[3],
                                                                    etay_1[4]))
        etay_mid = etay[5:-remainder]
        for k in np.arange(0,np.size(etay_mid), 5):
            records.append('      \t\t%.20e\t%.20e\t%.20e\t%.20e\t%.20e\n' % (etay_mid[k],
                                                                            etay_mid[k+1],
                                                                            etay_mid[k+2],
                                                                            etay_mid[k+3],
                                                                            etay_mid[k+4]))
        etay_last = etay[-remainder:]
        records.append('      \t\t%.20e' % etay_last[0])
        for e in etay_last[1:]:
            records.append('\t%.20e' %e)
        records.append('\n')
    elif n <= 5:
        etay = np.ones((n,1))*y[j-1]
        etay_1 = etay
        records.append('ETA=\t%i' % j)
        for e in etay_1:
            records.append('\t%.20e' % e)
        records.append('\n')

#Write to file
grd_fname = r'/Users/spmls/GitHub/pydelft/dev/test_data/write_rect_test.grd'
f = open(grd_fname, 'w')
for r in records:
    f.write(r)
f.close()

#print(records)

load grid and see if it works with depth file


In [10]:
grid = grd.grd()
grid.read_grd(grd_fname)
depth = dep.dep()
depth.read_dep(r'/Users/spmls/GitHub/pydelft/dev/test_data/test4.dep')

plt.pcolormesh(grid.x,grid.y,depth.depth)


Out[10]:
<matplotlib.collections.QuadMesh at 0x10cbf60f0>

test in pydelft now


In [15]:
m, n = 24, 14
Z = np.random.rand(m,n)
depth4 = dep.dep()
depth4.write_dep(Z, fname = r'/Users/spmls/GitHub/pydelft/dev/test_data/test_rectgrd.dep')
grid = grd.grd()
grid.write_rectgrd(fname = r'/Users/spmls/GitHub/pydelft/dev/test_data/test_rectgrd.grd', m = m, n = n, cellsize = 10,
                   x0 = 0, y0 = 0)

plt.pcolormesh(grid.x, grid.y, depth4.depth)


saved depth file: test_rectgrd.dep
Out[15]:
<matplotlib.collections.QuadMesh at 0x10ce29320>

In [12]: